home *** CD-ROM | disk | FTP | other *** search
/ C/C++ Users Group Library 1996 July / C-C++ Users Group Library July 1996.iso / vol_200 / 247_01 / bnpi.c < prev    next >
Encoding:
C/C++ Source or Header  |  1989-04-17  |  1.4 KB  |  67 lines
  1. /*
  2.  *   MIRACL calculate pi - by Gauss-Legendre method
  3.  *   bnpi.c
  4.  */
  5.   
  6. #include <stdio.h>
  7. #include "miracl.h"
  8.  
  9. /* Access global variables */
  10.  
  11. extern  int depth;     /* error tracing .. */
  12. extern  int trace[];   /* ... mechanism    */
  13. extern  int nib;       /* length of bigs   */
  14. extern  int lg2b;      /* bits in base     */
  15.  
  16. extern big  w11,w12,w13,w15; /* workspace variables */
  17.  
  18. void fpi(pi)
  19. flash pi;
  20. { /* Calculate pi using Guass-Legendre method */
  21.     int x,nits,op[5];
  22.     if (ERNUM) return;
  23.     depth++;
  24.     trace[depth]=53;
  25.     if (TRACER) track();
  26.     if (size(w15)!=0)
  27.     {
  28.         copy(w15,pi);
  29.         EXACT=FALSE;
  30.         depth--;
  31.         return;
  32.     }
  33.     fconv(1,2,pi);
  34.     froot(pi,2,pi);
  35.     fconv(1,1,w11);
  36.     fconv(1,4,w12);
  37.     x=1;
  38.     op[0]=0x6C;
  39.     op[1]=1;
  40.     op[4]=0;
  41.     nits=lg2b*nib/4;
  42.     while (x<nits)
  43.     {
  44.         copy(w11,w13);
  45.         op[2]=1;
  46.         op[3]=2;
  47.         flop(w11,pi,op,w11);
  48.         fmul(pi,w13,pi);
  49.         froot(pi,2,pi);
  50.         fsub(w11,w13,w13);
  51.         fmul(w13,w13,w13);
  52.         op[3]=1;
  53.         op[2]=(-x);
  54.         flop(w12,w13,op,w12);  /* w12 = w12 - x.w13 */
  55.         x*=2;
  56.     }
  57.     fadd(w11,pi,pi);
  58.     fmul(pi,pi,pi);
  59.     op[0]=0x48;
  60.     op[2]=0;
  61.     op[3]=4;
  62.     flop(pi,w12,op,pi);   /* pi = pi/(4.w12) */
  63.     copy(pi,w15);
  64.     depth--;
  65. }
  66.  
  67.